rm(list=ls())

##################################################################
### Campaign Finance Regulations and Public Policy -- Analysis ###
### Gilens, Patterson, and Haines 2020 ###########################
##################################################################

#############
# Libraries #
#############

library(panelView)
library(reshape2)
library(tidyverse)
library(usmap)

# Must install from devtools
# devtools::install_github('xuyiqing/gsynth')

library(gsynth)

#########################
# Set Working Directory #
#########################

#setwd('/Users/pattersos/Desktop/APSR Replication/Data/')

##################################
### Set Seed for Bootstrapping ###
##################################

set.seed(11121989)

###########################################
# Load in Dependent & Treatment Variables #
# Corporate Tax Rates & Ban Statuses ######
###########################################

ctax<-read.csv('state_corporate_tax.csv',stringsAsFactors = F)
names(ctax) <- c('state','year','ctax')
ctax<-subset(ctax,year<2017)

#############################
# Generate Treatment Status #
#############################

# Corp & Union Ban #
ctax$treat1<-ifelse(ctax$year>=2010 & ctax$state %in% c('AK','AZ','TX','OK','CO','WY','ND','SD','WI','MI','OH','PA','NC','NH','RI'), 1, 0)

# Corp Ban Only # 
ctax$treat2<-ifelse(ctax$year>=2010 &ctax$state %in% c('MT','MN','IA','KY','TN','WV','CT','MA'), 1, 0)

# At least Corp Ban # 
ctax$treat3<-ifelse(ctax$treat1==1 | ctax$treat2==1, 1, 0)

######################
# Load in Covariates #
######################

pop<-read.csv('state_population.csv')

### Economic Covariates ### 

# State GDP: https://www.bea.gov/data/gdp/gdp-state
sgdp<-read.csv('state_gdp.csv')
# State Unemployment: https://www.bls.gov/lau/#tables
unem<-read.csv('state_unemR.csv')
# State Deficits: https://slfdqs.taxpolicycenter.org/index.cfm
def<-read.csv('state_deficits.csv')

### Demand Covariates ### 

# State Firms: https://www.census.gov/data/tables/time-series/econ/susb/susb-historical.html
firms<-read.csv('state_firms.csv')
# State Unions: http://unionstats.gsu.edu/MonthlyLaborReviewArticle.htm
unions<-read.csv('state_unions.csv')
# State Individual Tax Rate: https://taxfoundation.org/state-individual-income-tax-rates/
indtax<-read.csv('state_individual_tax.csv')


### Political Covariates ### 

# Party Control: http://ukcpr.org/resources/national-welfare-data
party<-read.csv('state_party_control.csv')
# Party Vote Share: https://dataverse.harvard.edu/dataverse.xhtml?alias=cklarner
vote<-read.csv('state_party_vote.csv')
# Prez Vote Share: https://library.cqpress.com/elections/download-data.php
prez<-read.csv('state_presidential_vote.csv')

firms$state<-state.abb[match(firms$state,state.name)]
sgdp$state<-state.abb[match(sgdp$state,state.name)]
unem$state<-state.abb[match(unem$state,state.name)]
unions$state<-state.abb[match(unions$state,state.name)]

data<-merge(ctax,def,by=c('state','year'))
data<-merge(data,firms,by=c('state','year'))
data<-merge(data,indtax,by=c('state','year'))
data<-merge(data,party,by=c('state','year'))
data<-merge(data,prez,by=c('state','year'))
data<-merge(data,sgdp,by=c('state','year'))
data<-merge(data,unem,by=c('state','year'))
data<-merge(data,unions,by=c('state','year'))
data<-merge(data,vote,by=c('state','year'))
data<-merge(data,pop,by=c('state','year'))

#Which states don't have a corporate tax
#MI & OH, phase in / out structure 
#NV, SD, TX, WA, WY no corporate tax
#table(data$state[is.na(data$ctax)])

#################################################################
#### Alternative DVs ############################################
#################################################################
#################################################################
## Abortion Access (policy scale) ###############################
#################################################################

## Plicphys == must be preformed by a physician (1 = yes, 0 = no)
## Pgestlim == gestational limits (1 = yes, 0 = no)
## Pparbth == prohibit partial birth (1 = yes, 0 = no)
## Ppubfund == state funds beyond federal minimum (1 = yes, 0 = no)
## Pprivins == private insurance restrictions (1 = yes, 0 = no)
## Pwait == mandate waiting period (1 = yes, 0 = no)
## Pminors == Pnotif + Pcons (parental notification / consent)

abort<-read.csv('state_abortion_regulations.csv',stringsAsFactors = F)
abort$State<-state.abb[match(abort$State,state.name)]


# use 30-item index of gun laws from Sorens dataset
# scaled to run from 0 to 100 with high scores=liberal (stricter) gun laws
# btrain    Brifle    Bmult     bssp      Blicens   Bballist  Rbsi      RBpreoc   RBbmg     Rbsbs
# Basslt    Bminage   Bdealer   botheft   Bregis    Bretent   RBopen    RBprecc   Rbssp     RBjourn
# Bmags     Bwait     bstheft   Bpriv     Bdesign   Rbccpir   RBplaces  Rbdtr     Rbsbr     RBlocks

guns<-read.csv('state_gun_laws.csv',stringsAsFactors = F)
guns$state<-state.abb[match(guns$state,state.name)]

#################################################################
## Eminent Domain and Utility Regulations 
#################################################################
# use Marty's new eminent domain that excludes the consitutional item
# high scores reflect more restrictive, anti-corporate laws (state can't take land for development)
# emdom3 = (freform + fprivate + fblight)/3
#################################################################

land<-read.csv('state_emdom_laws.csv',stringsAsFactors = F)
land$State<-state.abb[match(land$State,state.name)]

## freform == enacted eminent domain reform, stricter limits on "public use" 
## fprivate == prohibits private taking for private use (1 = all, .5 = some, 0 = none)
## fblight == enacted stricter defintion of 'blight' (1 = yes, .5 = somewhat, 0 = no)
## fconst == constitutoin enshrines additional restrictions (1 = yes, 0 = no, .5 = some)
## findex == eminent domain index ((freform+fprivate+fblight)*1+(.5*fconst))
## frenew == mandatory renewable portfolio standard 
## frpsexc == renewable portfolio standard exceptions (1.5 = hydro, .5 gas or nuclear, 1 = otherwise)
## frps == renewable porfolio index

#################################################################
## Civil Liabilities / Plaintiff Friendliness of Legal System
#############################################################
## use index of 3 tort items; index created by Marty in Stata
## average of rescaled opuncap, opunburd, ojsl
## scored 0-1 with high score being anti-corporate & plaintiff-friendly
#############################################################

civil<-read.csv('state_tort_laws.csv',stringsAsFactors = F)
civil$State<-state.abb[match(civil$State,state.name)]
names(civil)[1]<-'civil'

## oliamba == three year smoothed tort system score
## opartsc == partisan supreme court (1 = yes, 0 = no, .1 = retention)
## opartrial == partisan trial courts 
## olawyers == lawyer concentraion
## opuncap == punative / noneconomic damages cap (2 = all prohibited, 1 = yes, 0 = no)
## opunburd == evidence standard for pun burd (2 = reasonable d, 1 = clear and con, 0 = prepon, 3 = no putative)
## ojsl == joint and several liability abolished (1 = yes, .5 = some, 0 = no)
## ocourts == first principal component, plaintiff-friendly civil liability index

#############################################################

abort<-abort[c('State','Year','abortion')]
guns<-guns[c('state','year','guns30')]
civil<-civil[c('State','Year','civil')]
land<-land[c('State','Year','emdom3')]

names(abort)<-tolower(names(abort))
names(civil)<-tolower(names(civil))
names(land)<-tolower(names(land))

data$firmsPC<-data$firms500/(data$population/100000)
data$ctax<-data$ctax*100

data2<-merge(data,guns,by=c('state','year'),all.x=T)
data2<-merge(data2,abort,by=c('state','year'),all.x=T)
data2<-merge(data2,civil,by=c('state','year'),all.x=T)
data2<-merge(data2,land,by=c('state','year'),all.x=T)

setdiff(names(data2),names(data))

rm(list=setdiff(ls(), c("data",'data2')))

cu<-na.omit(data)

## Central Model
## Corporate Taxes 
## Interactive Fixed Effects

## recode ctax to 0-100 range

min(cu$ctax)
max(cu$ctax)
cu$ctax100 <- ((cu$ctax-3.4)/8.6)*100
min(cu$ctax100)
max(cu$ctax100)

out_taxA<-gsynth(ctax100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,
                 data=cu, EM=T, index=c('state','year'),force='two-way', 
                 CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,
                 parallel=F,cores=1, seed = 100)
out_taxA;plot(out_taxA, main = "Higher Corporate Tax Rates - All Affected States")


## IFE for ban only corporate spending states (ie, include only no ban and corp ban states)

both<-c('AK','AZ','TX','OK','CO','WY','ND','SD','WI','MI','OH','PA','NC','NH','RI')
corp_only<-c('MT','MN','IA','KY','TN','WV','CT','MA')
treated<-c(both,corp_only)
untreated<-setdiff(state.abb,treated)

cu2 <- cu %>% filter(state %in% untreated | state %in% corp_only)
cu2 <- na.omit(cu2)

out_taxB<-gsynth(ctax100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,
               data=cu2, EM=T, index=c('state','year'),force='two-way', 
               CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,
               parallel=F,cores=1, seed = 100)
out_taxB;plot(out_taxB, main = "Higher Corporate Tax Rates - Corp Ban States Only")


#######################
### Abortion Models ###
#######################

abortion<- data2 %>% select(govD,abortion,state,year,treat3,Dgdp,unemR,PCdeficit16,firmsPC,union,indtax,rep_share,Rprez,party)
abortion<-na.omit(abortion)

## Rescale to 0-100 and so higher scores are more liberal (to parallel higher taxes) 
min(abortion$abortion)
max(abortion$abortion)
abortion$abortion100 <- ((abortion$abortion*100)*-1)+100

# All Treated 
out_abortionA<-gsynth(abortion100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,
                      data=abortion, EM=T, index=c('state','year'), force='two-way', 
                      CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,
                      parallel=F,cores=1, seed = 100)
out_abortionA;plot(out_abortionA, main = "Less Restrictive Abortion Laws - All Affected States")

# Just Corporate Bans
ab2 <- abortion %>% filter(state %in% untreated | state %in% corp_only)
ab2 <- na.omit(ab2)


out_abortionB<-gsynth(abortion100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,
                      data=ab2, EM=T, index=c('state','year'), force='two-way', 
                      CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,
                      parallel=F,cores=1, seed = 100)
out_abortionB;plot(out_abortionB, main = "Less Restrictive Abortion Laws - Corp Ban States Only")

##########################
### Gun Control Models ###
##########################

guns <- data2 %>% select(govD,guns30,state,year,treat3,Dgdp,unemR,PCdeficit16,firmsPC,union,indtax,rep_share,Rprez,party)
guns<-na.omit(guns)

# All Treated 
out_gunsA<-gsynth(guns30~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,
                      data=guns, EM=T, index=c('state','year'), force='two-way', 
                      CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,
                      parallel=F,cores=1, seed = 100)
out_gunsA;plot(out_gunsA, main = "Strong Gun Control Laws - All Affected States")

# Just Corporate Bans
guns2 <- guns %>% filter(state %in% untreated | state %in% corp_only)
guns2 <- na.omit(guns2)


out_gunsB<-gsynth(guns30~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,
                      data=guns2, EM=T, index=c('state','year'), force='two-way', 
                      CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,
                      parallel=F,cores=1, seed = 100)
out_gunsB;plot(out_gunsB, main = "Strong Gun Control Laws - Corp Ban States Only")

##############################
### Emminent Domain Models ###
##############################

emdom <- data2 %>% select(govD,emdom3,state,year,treat3,Dgdp,unemR,PCdeficit16,firmsPC,union,indtax,rep_share,Rprez,party)
emdom<-na.omit(emdom)

## rescale to 0-100
emdom$emdom3100 <- emdom$emdom3*100

# All Treated 
out_emdomA<-gsynth(emdom3100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,
                  data=emdom, EM=T, index=c('state','year'), force='two-way', 
                  CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,
                  parallel=F,cores=1, seed = 100)
out_emdomA;plot(out_emdomA, main = "Anti-Corporate Eminent Domain Laws - All Affected States")

# Just Corporate Bans
emdom2 <- emdom %>% filter(state %in% untreated | state %in% corp_only)
emdom2 <- na.omit(emdom2)


out_emdomB<-gsynth(emdom3100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,
                  data=emdom2, EM=T, index=c('state','year'), force='two-way', 
                  CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,
                  parallel=F,cores=1, seed = 100)
out_emdomB;plot(out_emdomB, main = "Anti-Corporate Eminent Domain Laws - Corp Ban States Only")

#######################
### Tort Law Models ###
#######################

civil <- data2 %>% select(govD,civil,state,year,treat3,Dgdp,unemR,PCdeficit16,firmsPC,union,indtax,rep_share,Rprez,party)
civil<-na.omit(civil)

## rescale to 0-100
civil$civil100 <- civil$civil*100

# All Treated 
out_civilA<-gsynth(civil100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,
                   data=civil, EM=T, index=c('state','year'), force='two-way', 
                   CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,
                   parallel=F,cores=1, seed = 100)
out_civilA;plot(out_civilA, main = "Plaintiff-Friendly Tort Laws - All Affected States")

# Just Corporate Bans
civil2 <- civil %>% filter(state %in% untreated | state %in% corp_only)
civil2 <- na.omit(civil2)


out_civilB<-gsynth(civil100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,
                   data=civil2, EM=T, index=c('state','year'), force='two-way', 
                   CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,
                   parallel=F,cores=1, seed = 100)
out_civilB;plot(out_civilB, main = "Plaintiff-Friendly Tort Laws - Corp Ban States Only")

##########
# Review #
##########

out_taxA
fig_taxA<-plot(out_taxA)
data_figA <- fig_taxA$data
taxfigA <- data_figA %>% ggplot(aes(x=time,y=ATT)) + geom_hline(yintercept = 0,lwd=1.25,col='grey80') + geom_vline(xintercept = 0,lwd=1.25,col='grey80') + 
  geom_line(lwd=1.5) + geom_ribbon(aes(x=time,ymin=CI.lower,ymax=CI.upper),na.rm=F,orientation = NA,outline.type = 'both',alpha=.2, flipped_aes=F) + 
  theme_bw() + xlab('Years Since Citizens United') + ylab("Estimated Effect") + 
  ggtitle('Higher Corporate Tax Rates',subtitle = 'All Affected States')

out_taxB
fig_taxB<-plot(out_taxB)
data_figB <- fig_taxB$data
taxfigB <- data_figB %>% ggplot(aes(x=time,y=ATT)) + geom_hline(yintercept = 0,lwd=1.25,col='grey80') + geom_vline(xintercept = 0,lwd=1.25,col='grey80') + 
  geom_line(lwd=1.5) + geom_ribbon(aes(x=time,ymin=CI.lower,ymax=CI.upper),na.rm=F,orientation = NA,outline.type = 'both',alpha=.2, flipped_aes=F) + 
  theme_bw() + xlab('Years Since Citizens United') + ylab("Estimated Effect") + 
  ggtitle('Higher Corporate Tax Rates',subtitle = 'Corporate Ban Only')


out_abortionA
fig_abortionA<-plot(out_abortionA)
data_figA <- fig_abortionA$data
abortionfigA <- data_figA %>% ggplot(aes(x=time,y=ATT)) + geom_hline(yintercept = 0,lwd=1.25,col='grey80') + geom_vline(xintercept = 0,lwd=1.25,col='grey80') + 
  geom_line(lwd=1.5) + geom_ribbon(aes(x=time,ymin=CI.lower,ymax=CI.upper),na.rm=F,orientation = NA,outline.type = 'both',alpha=.2, flipped_aes=F) + 
  theme_bw() + xlab('Years Since Citizens United') + ylab("Estimated Effect") + 
  ggtitle('Less Restrictive Abortion Laws',subtitle = 'All Affected States')

out_abortionB
fig_abortionB<-plot(out_abortionB)
data_figB <- fig_abortionB$data
abortionfigB <- data_figB %>% ggplot(aes(x=time,y=ATT)) + geom_hline(yintercept = 0,lwd=1.25,col='grey80') + geom_vline(xintercept = 0,lwd=1.25,col='grey80') + 
  geom_line(lwd=1.5) + geom_ribbon(aes(x=time,ymin=CI.lower,ymax=CI.upper),na.rm=F,orientation = NA,outline.type = 'both',alpha=.2, flipped_aes=F) + 
  theme_bw() + xlab('Years Since Citizens United') + ylab("Estimated Effect") + 
  ggtitle('Less Restrictive Abortion Laws',subtitle = 'Corporate Ban Only')

out_gunsA
fig_gunsA<-plot(out_gunsA)
data_figA <- fig_gunsA$data
gunsfigA <- data_figA %>% ggplot(aes(x=time,y=ATT)) + geom_hline(yintercept = 0,lwd=1.25,col='grey80') + geom_vline(xintercept = 0,lwd=1.25,col='grey80') + 
  geom_line(lwd=1.5) + geom_ribbon(aes(x=time,ymin=CI.lower,ymax=CI.upper),na.rm=F,orientation = NA,outline.type = 'both',alpha=.2, flipped_aes=F) + 
  theme_bw() + xlab('Years Since Citizens United') + ylab("Estimated Effect") + 
  ggtitle('Strong Gun Control Laws',subtitle = 'All Affected States')

out_gunsB
fig_gunsB<-plot(out_gunsB)
data_figB <- fig_gunsB$data
gunsfigB <- data_figB %>% ggplot(aes(x=time,y=ATT)) + geom_hline(yintercept = 0,lwd=1.25,col='grey80') + geom_vline(xintercept = 0,lwd=1.25,col='grey80') + 
  geom_line(lwd=1.5) + geom_ribbon(aes(x=time,ymin=CI.lower,ymax=CI.upper),na.rm=F,orientation = NA,outline.type = 'both',alpha=.2, flipped_aes=F) + 
  theme_bw() + xlab('Years Since Citizens United') + ylab("Estimated Effect") + 
  ggtitle('Strong Gun Control Laws',subtitle = 'Corporate Ban Only')

out_emdomA
fig_emdomA<-plot(out_emdomA)
data_figA <- fig_emdomA$data
emdomfigA <- data_figA %>% ggplot(aes(x=time,y=ATT)) + geom_hline(yintercept = 0,lwd=1.25,col='grey80') + geom_vline(xintercept = 0,lwd=1.25,col='grey80') + 
  geom_line(lwd=1.5) + geom_ribbon(aes(x=time,ymin=CI.lower,ymax=CI.upper),na.rm=F,orientation = NA,outline.type = 'both',alpha=.2, flipped_aes=F) + 
  theme_bw() + xlab('Years Since Citizens United') + ylab("Estimated Effect") + 
  ggtitle('Anti-Corporate Eminent Domain Laws',subtitle = 'All Affected States')

out_emdomB
fig_emdomB<-plot(out_emdomB)
data_figB <- fig_emdomB$data
emdomfigB <- data_figB %>% ggplot(aes(x=time,y=ATT)) + geom_hline(yintercept = 0,lwd=1.25,col='grey80') + geom_vline(xintercept = 0,lwd=1.25,col='grey80') + 
  geom_line(lwd=1.5) + geom_ribbon(aes(x=time,ymin=CI.lower,ymax=CI.upper),na.rm=F,orientation = NA,outline.type = 'both',alpha=.2, flipped_aes=F) + 
  theme_bw() + xlab('Years Since Citizens United') + ylab("Estimated Effect") + 
  ggtitle('Anti-Corporate Eminent Domain Laws',subtitle = 'Corporate Ban Only')

out_civilA
fig_civilA<-plot(out_civilA)
data_figA <- fig_civilA$data
civilfigA <- data_figA %>% ggplot(aes(x=time,y=ATT)) + geom_hline(yintercept = 0,lwd=1.25,col='grey80') + geom_vline(xintercept = 0,lwd=1.25,col='grey80') + 
  geom_line(lwd=1.5) + geom_ribbon(aes(x=time,ymin=CI.lower,ymax=CI.upper),na.rm=F,orientation = NA,outline.type = 'both',alpha=.2, flipped_aes=F) + 
  theme_bw() + xlab('Years Since Citizens United') + ylab("Estimated Effect") + 
  ggtitle('Plaintiff-Friendly Tort Laws',subtitle = 'All Affected States')

out_civilB
fig_civilB<-plot(out_civilB)
data_figB <- fig_civilB$data
civilfigB <- data_figB %>% ggplot(aes(x=time,y=ATT)) + geom_hline(yintercept = 0,lwd=1.25,col='grey80') + geom_vline(xintercept = 0,lwd=1.25,col='grey80') + 
  geom_line(lwd=1.5) + geom_ribbon(aes(x=time,ymin=CI.lower,ymax=CI.upper),na.rm=F,orientation = NA,outline.type = 'both',alpha=.2, flipped_aes=F) + 
  theme_bw() + xlab('Years Since Citizens United') + ylab("Estimated Effect") + 
  ggtitle('Plaintiff-Friendly Tort Laws',subtitle = 'Corporate Ban Only')


#pdf('/users/Shawn/desktop/Fig1A.pdf',height=5,width=10)
gridExtra::grid.arrange(taxfigA,taxfigB,ncol=2)
#dev.off()

#pdf('/users/Shawn/desktop/Fig1B.pdf',height=10,width=10)
gridExtra::grid.arrange(abortionfigA,abortionfigB,gunsfigA,gunsfigB,ncol=2)
#dev.off()

#pdf('/users/Shawn/desktop/Fig1C.pdf',height=10,width=10)
gridExtra::grid.arrange(emdomfigA,emdomfigB,civilfigA,civilfigB,ncol=2)
#dev.off()

##################
### Appendix C ###
##################

out_taxA;out_taxB
out_civilA;out_civilB
out_emdomA;out_emdomB
out_abortionA;out_abortionB
out_gunsA;out_gunsB 

##################
### Appendix D ###
##################

out_taxC <-gsynth(ctax100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+rep_share+Rprez,data=cu, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)
out_taxD <-gsynth(ctax100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+rep_share+Rprez,data=cu2, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)

out_civilC <- gsynth(civil100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+rep_share+Rprez,data=civil, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)
out_civilD <- gsynth(civil100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+rep_share+Rprez,data=civil2, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)

out_emdomC <- gsynth(emdom3100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+rep_share+Rprez,data=emdom, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)
out_emdomD <- gsynth(emdom3100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+rep_share+Rprez,data=emdom2, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)

out_abortionC <- gsynth(abortion100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+rep_share+Rprez,data=abortion, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)
out_abortionD <- gsynth(abortion100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+rep_share+Rprez,data=ab2, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)

out_gunsC <- gsynth(guns30~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+rep_share+Rprez,data=guns, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)
out_gunsD <- gsynth(guns30~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+rep_share+Rprez,data=guns2, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)

out_taxC;out_taxD
out_civilC;out_civilD
out_emdomC;out_emdomD
out_abortionC;out_abortionD
out_gunsC;out_gunsD 

##################
### Appendix E ###
##################

out_taxE <-gsynth(ctax100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,data=cu, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, estimator = 'mc', inference='nonparametric',nboots=1000,parallel=F,cores=1, seed = 100)
out_taxF <-gsynth(ctax100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,data=cu2, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, estimator = 'mc', inference='nonparametric',nboots=1000,parallel=F,cores=1, seed = 100)

out_civilE <- gsynth(civil100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,data=civil, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, estimator = 'mc', inference='nonparametric',nboots=1000,parallel=F,cores=1, seed = 100)
out_civilF <- gsynth(civil100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,data=civil2, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, estimator = 'mc', inference='nonparametric',nboots=1000,parallel=F,cores=1, seed = 100)

out_emdomE <- gsynth(emdom3100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,data=emdom, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, estimator = 'mc', inference='nonparametric',nboots=1000,parallel=F,cores=1, seed = 100)
out_emdomF <- gsynth(emdom3100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,data=emdom2, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, estimator = 'mc', inference='nonparametric',nboots=1000,parallel=F,cores=1, seed = 100)

out_abortionE <- gsynth(abortion100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,data=abortion, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, estimator = 'mc', inference='nonparametric',nboots=1000,parallel=F,cores=1, seed = 100)
out_abortionF <- gsynth(abortion100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,data=ab2, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, estimator = 'mc', inference='nonparametric',nboots=1000,parallel=F,cores=1, seed = 100)

out_gunsE <- gsynth(guns30~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,data=guns, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, estimator = 'mc', inference='nonparametric',nboots=1000,parallel=F,cores=1, seed = 100)
out_gunsF <- gsynth(guns30~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,data=guns2, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, estimator = 'mc', inference='nonparametric',nboots=1000,parallel=F,cores=1, seed = 100)


out_taxE;out_taxF
out_civilE;out_civilF
out_emdomE;out_emdomF
out_abortionE;out_abortionF
out_gunsE;out_gunsF

##################
### Appendix F ###
##################

treated<-c('AK','AZ','CO','CT','IA','KY','MA','MN','MT','NC','ND','NH','OK','PA','RI','TN','WI','WV')


for(i in 1:length(treated)){
  
  cuP <- cu %>% filter(state!=treated[i])
  
  out_tax_out<-gsynth(ctax100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,
                   data=cuP, EM=T, index=c('state','year'),force='two-way', 
                   CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,
                   parallel=F,cores=1)
  
  print(out_tax_out)
  assign(paste(treated[i]),out_tax_out)
  print(paste(treated[i]))
}

##################
### Appendix G ###
##################

out_taxG <-gsynth(ctax100~treat3,data=cu, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)

out_civilG <- gsynth(civil100~treat3,data=civil, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)
out_emdomG <- gsynth(emdom3100~treat3,data=emdom, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)
out_abortionG <- gsynth(abortion100~treat3,data=abortion, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)
out_gunsG <- gsynth(guns30~treat3,data=guns, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)

nocorp <- setdiff(unique(abortion$state),unique(cu$state))

ab3 <- abortion %>% filter(!(state %in% nocorp))
guns3 <- guns %>% filter(!(state %in% nocorp))
emdom3 <- emdom %>% filter(!(state %in% nocorp))
civil3 <- civil %>% filter(!(state %in% nocorp))

out_civilH <- gsynth(civil100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,data=civil3, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T,  inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)
out_emdomH <- gsynth(emdom3100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,data=emdom3, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T,  inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)
out_abortionH <- gsynth(abortion100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,data=ab3, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T,  inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)
out_gunsH <- gsynth(guns30~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+party+Rprez,data=guns3, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)


out_taxG
out_civilG;out_civilH
out_emdomG;out_emdomH
out_abortionG;out_abortionH
out_gunsG;out_gunsH  

##################
### Appendix H ###
##################

out_taxI <-gsynth(ctax100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+govD+party+Rprez,data=cu, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T,  inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)

out_civilI <- gsynth(civil100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+govD+party+Rprez,data=civil, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T,  inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)

out_emdomI <- gsynth(emdom3100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+govD+party+Rprez,data=emdom, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T,  inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)

out_abortionI <- gsynth(abortion100~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+govD+party+Rprez,data=abortion, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T,  inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)

out_gunsI <- gsynth(guns30~treat3+Dgdp+unemR+PCdeficit16+firmsPC+union+indtax+govD+party+Rprez,data=guns, EM=T, index=c('state','year'),force='two-way', CV=T,r=c(0,5), se=T, inference='parametric',nboots=1000,parallel=F,cores=1, seed = 100)

out_taxI
out_civilI
out_emdomI
out_abortionI
out_gunsI

##################
### Appendix I ###
##################


both<-c('AK','AZ','TX','OK','CO','WY','ND','SD','WI','MI','OH','PA','NC','NH','RI')
corp_only<-c('MT','MN','IA','KY','TN','WV','CT','MA')
statepop$ban <- as.factor(ifelse(statepop$abbr %in% corp_only,1,ifelse(statepop$abbr %in% both,2,0)))

map <- plot_usmap(data = statepop, values = "ban", color = "black") +
  scale_fill_manual('I.E. Bans',values = c("white", "grey70", "grey30"), 
                    label=c('No I.E. Ban','Corporate Only','Corporate & Union')) + theme(legend.position = c(.6,0))

#pdf('/users/shawn/Dropbox/State policy/writing/APSR revisions/Map.pdf',height=6,width=8)
print(map)
#dev.off()




g_sum<-function(x,y){
  temp1<-x$est.avg[,c(1,2,5)]
  temp2<-x$est.att[,c(1,2,5,6)]
  
  temp3<-y$est.avg[,c(1,2,5)]
  temp4<-y$est.att[,c(1,2,5,6)]
  
  mat <- matrix(NA,nrow=nrow(temp2)+1,ncol=9)

  mat[2:nrow(mat),2:5]<-temp2
  mat[1,2:4]<-temp1 
  
  mat[2:nrow(mat),6:9]<-temp4
  mat[1,6:8]<-temp3 
  
  mat[,1] <- c('Avg.',rev(seq(2016,by=-1,length=nrow(temp2))))
  mat <- as.data.frame(mat)
  
  mat[1,5]<-mat[nrow(mat),5]
  mat[1,9]<-mat[nrow(mat),9]
  
  names(mat)<-c('Year','ATT','S.E.','p-value','N','ATT2','S.E.2','p-value2','N2')
  return(mat)

}


write_csv(g_sum(out_abortionA,out_abortionB),'/Users/pattersos/Desktop/APSR Replication/Appendix Tables/abortion_AppC.csv')
write_csv(g_sum(out_abortionC,out_abortionD),'/Users/pattersos/Desktop/APSR Replication/Appendix Tables/abortion_AppD.csv')
write_csv(g_sum(out_abortionE,out_abortionF),'/Users/pattersos/Desktop/APSR Replication/Appendix Tables/abortion_AppE.csv')

write_csv(g_sum(out_civilA,out_civilB),'/Users/pattersos/Desktop/APSR Replication/Appendix Tables/civil_AppC.csv')
write_csv(g_sum(out_civilC,out_civilD),'/Users/pattersos/Desktop/APSR Replication/Appendix Tables/civil_AppD.csv')
write_csv(g_sum(out_civilE,out_civilF),'/Users/pattersos/Desktop/APSR Replication/Appendix Tables/civil_AppE.csv')

write_csv(g_sum(out_emdomA,out_emdomB),'/Users/pattersos/Desktop/APSR Replication/Appendix Tables/emdom_AppC.csv')
write_csv(g_sum(out_emdomC,out_emdomD),'/Users/pattersos/Desktop/APSR Replication/Appendix Tables/emdom_AppD.csv')
write_csv(g_sum(out_emdomE,out_emdomF),'/Users/pattersos/Desktop/APSR Replication/Appendix Tables/emdom_AppE.csv')

write_csv(g_sum(out_gunsA,out_gunsB),'/Users/pattersos/Desktop/APSR Replication/Appendix Tables/guns_AppC.csv')
write_csv(g_sum(out_gunsC,out_gunsD),'/Users/pattersos/Desktop/APSR Replication/Appendix Tables/guns_AppD.csv')
write_csv(g_sum(out_gunsE,out_gunsF),'/Users/pattersos/Desktop/APSR Replication/Appendix Tables/guns_AppE.csv')

write_csv(g_sum(out_taxA,out_taxB),'/Users/pattersos/Desktop/APSR Replication/Appendix Tables/tax_AppC.csv')
write_csv(g_sum(out_taxC,out_taxD),'/Users/pattersos/Desktop/APSR Replication/Appendix Tables/tax_AppD.csv')
write_csv(g_sum(out_taxE,out_taxF),'/Users/pattersos/Desktop/APSR Replication/Appendix Tables/tax_AppE.csv')




